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ABSTRACT: 

This paper presents the design and implementation 
of the tetrolets based system for automatic skeletal Bone 
Age Assessment (BAA). The system works according to 
the renowned Tanner and Whitehouse (TW2) method, 
based on the carpal and phalangeal Region of Interest 
(ROI). The system ensures accurate and robust BAA for 
the age range 0-10 years for both girls and boys. Given 
a left hand-wrist radiograph as input, the system 
estimates the bone age by deploying novel techniques 
for segmentation, feature extraction, feature selection 
and classification. Tetrolets are used in combination 
with Particle Swarm Optimization (PSO) for 
segmentation. From the segmented wrist bones, the 
carpal and phalangeal ROI are identified and are used in 
morphological feature extraction. PCA is employed as a 
feature selection tool to reduce the size of the feature 
vector. The selected features are fed in to an ID3 
decision tree classifier, which outputs the class to which 
the radiograph is categorized, which is mapped onto the 
final bone age. The system was evaluated on a set of 
100 radiographs (50 for girls and 50 for boys), and the 
results are discussed. The performance of system was 
evaluated with the help of radiologist expert diagnoses. 
The system is very reliable with minimum human 
intervention, yielding excellent results. 

Keywords: Bone Age Assessment (BAA), TW2, 
radiograph. Particle Swarm Optimization (PSO), 
Tetrolets, ID3, Classification. 

1. Introduction: 

The chronological situations of humans are 
described by certain indices such as height, dental age, 
and bone maturity. Of these, bone age measurement 
plays a significant role because of its reliability and 
practicability in diagnosing hereditary diseases and 
growth disorders. Bone age assessment using a hand 
radiograph is an important clinical tool in the area of 
pediatrics, especially in relation to endocrinological 
problems and growth disorders. A single reading of 
skeletal age informs the clinician of the relative maturity 
of a patient at a particular time in his or her life and 
integrated with other clinical finding, separates the 


normal from the relatively advanced or retarded [1]. 
The bone age of children is apparently influenced by 
gender, race, nutrition status, living environments and 
social resources, etc. Based on a radiological 
examination of skeletal development of the left-hand 
wrist, bone age is assessed and compared with the 
chronological age. A discrepancy between these two 
values indicates abnormalities in skeletal development. 
The procedure is often used in the management and 
diagnosis of endocrine disorders and also serves as an 
indication of the therapeutic effect of treatment. It 
indicates whether the growth of a patient is accelerating 
or decreasing, based on which the patient can be treated 
with growth hormones. BAA is universally used due to 
its simplicity, minimal radiation exposure, and the 
availability of multiple ossification centers for 
evaluation of maturity. 


2. Background of BAA: 

The main clinical methods for skeletal bone age 
estimation are the Greulich & Pyle (GP) method and the 
Tanner & Whitehouse (TW) method. GP is an atlas 
matching method while TW is a score assigning method 
[2]. GP method is faster and easier to use than the TW 
method. Bull et. al. performed a large scale comparison 
of the GP and TW method and concluded that TW 
method is the more reproducible of the two and 
potentially more accurate 13]. 
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Fig. 1 . Bones of hand and wrist for BAA 
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In GP method, a left-hand wrist radiograph is 
compared with a series of radiographs grouped in the 
atlas according to age and sex. The atlas pattern which 
superficially appears to resemble the clinical image is 
selected. Since each atlas pattern is assigned to a certain 
year of age, the selection assesses the bone age. The 
disadvantage of this method is the subjective nature of 
the analysis performed by various observers with 
different levels of training. The reason for high 
discrepancies in atlas matching method is due to a 
general comparison of the radiograph to the atlas 
pattern. By a more detailed comparison of individual 
bones, ambiguous results may be obtained. 

TW method uses a detailed analysis of each 
individual bone (shown in Fig. 1), assigning it to one of 
eight classes reflecting its developmental stage. This 
leads to the description of each bone in terms of scores. 
The sum of all scores assesses the bone age. This 
method yields the most reliable results. The high 
complexity of the TW method is the main reason for its 
less intensive use and what makes it worthwhile to 
automate. The original Tanner- Whitehouse method 
(TW1, 1962) was presented by Tanner, Whitehouse and 
Healy. The fundamental advantage of TW 1 was its solid 
and formal mathematical soul. Based on stages, scores 
were assigned and later on added to obtain the final 
skeletal age. TW2 was a revision of TW1, especially in 
relation to the scores associated to each stage and also 
the difference between both sexes. The TW2 method 
does not use a scale based on the age, rather it is based 
on a set of bone’s standard maturity for each age 
population. In detail, in the TW2 method twenty regions 
of interest (ROIs), located in the main bones are 
considered for the bone age evaluation. Each ROI is 
divided into three parts: Epiphysis, Metaphysis and 
Diaphysis; it is possible to identify these different 
ossification centers in the phalanx proximity. The 
development of each ROI is divided into discrete stages, 
as shown in Fig. 2, and each stage is given a letter 
(A,B,C,D,...I), reflecting the development stage as: 

■ Stage A - absent 

■ Stage B - single deposit of calcium 

■ Stage C - center is distinct in appearance 


■ Stage D - maximum diameter is half or more the 
width of metaphysis 

■ Stage E - border of the epiphysis is concave 

■ Stage F - epiphysis is as wide as metaphysis 

■ Stage G - epiphysis caps the metaphysis 

■ Stage H - fusion of epiphysis and metaphysis has 
begun 

■ Stage 1 - epiphyseal fusion completed. 

By adding the scores of all ROIs, an overall maturity 
score is obtained. This score is correlated with the bone 
age differently for males and females [4]. For TW2 
method, these score systems have been developed: 

■ TW2 20 bones: characterized by 20 bones including 
the bones of the first, third and fifth finger and the 
carpal bones. 

■ RUS: considers the same bones of the TW2 method 
except the carpal bones. 

■ CARPAL: considers only the carpal bones. 

A number of algorithms for automated skeletal bone age 
assessment exist in the literature. 

3. Survey of Literature: 

In early 1980s, Pal and King proposed the theory of 
fuzzy sets and applied it for edge detection algorithm of 
X-ray images 15]. Kwabwe et. al. later in 1986, 
proposed certain algorithms to recognize the bones in an 
X-ray image of the hand and wrist [6]. They used a 
shape description technique based on linear 
measurements from a polygonal approximation of the 
bones. A fuzzy classifier for syntactic recognition of 
different stages of maturity of bones from X-rays of 
hand and wrist using fuzzy grammar and fuzzy 
primitives was developed by Pathak and Pal [7]. It 
comprised of a hierarchical three-stage syntactic 
recognition algorithm, which made use of six-tuple 
fuzzy and seven-tuple fractionally fuzzy grammars to 
identify the different stages of maturity of bones from 
X-rays. Michael and Nelson [8] developed a model- 
based system for automatic segmentation of bones from 
digital hand radiographs named as HANDX, in 1989. 
This computer vision system, offered a solution to 
automatically find, isolate and measure bones from 
digital X-rays. In 1991, Pietka et. al. described a method 
19] based on independent analysis of the phalangeal 
regions. Phalangeal analysis was performed in several 
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stages by measuring the lengths of the distal, middle and 
proximal phalanx. These measurements were converted 
into skeletal age by using the standard phalangeal length 
table proposed by Gam et.al [1()|. These single bone age 
estimates were then averaged to assess the global 
phalangeal age of the patient. 

Tanner and Gibbons introduced the Computer- 
Assisted Skeletal Age Scores (CASAS) system in 1992 
[111- This was based on nine prototype images for each 
bone, representing the nine stages of maturity. Thus, a 
stage was defined by an image template. Two or three 
most similar templates for the radiographs were 
identified. The system then automatically computed a 
measure of correlation to each template and a fractional 
stage. The correlation to the template was a measure of 
similarity. In 1993, Pietka et. al. performed phalangeal 
and carpal bone analysis using standard and dynamic 
thresholding methods to assess skeletal age [12]. Cheng 
et. al. [13] proposed the methods to extract a region of 
interest (ROI) for texture analysis in 1994, with 
particular attention to patients with 
hyperparathyroidism. The techniques included 
multiresolution sensing, automatic adaptive 
thresholding, detection of orientation angle, and 
projection taken perpendicular to the line of least second 
moment. In the same year, Drayer and Cox [14] 
designed a computer aided system to estimate bone age 
based on Fourier analysis on radiographs to produce 
TW2 standards for radius, ulna and short finger bones. 
It employed template matching of each bone to the 
scanned image of the radiograph. In 1996, Al-Taani et. 
al. classified the bones of the hand -wrist images into 
pediatric stages of maturity using Point Distribution 
Models (PDM) [15]. Wastl and Dickhaus proposed a 
pattern recognition based BAA approach, in the same 
year [16]. The approach consisted of four major steps: 
digitization of the hand radiograph, segmentation of 
ROI, prototype matching and BAA. In 1999, Bull et. al. 
made a remarkable comparison of GP and TW2 
methods [3] and concluded the following. The GP 
method involves a complex comparison of all of the 
bones in the hand and wrist against reference “normal” 
radiographs of different ages. Although this approach is 
considerably faster than the original, it may be less 
accurate. The TW2 method relies on the systematic 
evaluation of the maturity of all the bones in the hand 
and wrist. The measured intra-observer variation was 
greater for the GP method than for the TW2 method. 
This accounts for much of the discrepancy between the 
two methods. They also concluded that the GP and TW2 
methods produce different values for bone age, which 
are significant in clinical practice. They have also 
shown that the TW2 method is more reproducible than 


the GP method. They finally suggested TW2 method to 
be preferably used as the only one BAA method when 
performing serial measurements of a patient. Mahmoodi 
et. al. (1997) used Knowledge-based Active Shape 
Models (ASM) in an automated vision system to assess 
the bone age [17]. Pietka et. al. conducted a computer 
assisted BAA procedure [18] by extracting and using 
the epiphyseal/ metaphyseal ROI (EMROI), in 2001. 
From each phalanx 3 EMROIs were extracted which 
include: metaphysis, epiphysis and diaphysis of the 
distal and middle phalanges and for the proximal 
phalanges it includes metaphysis, epiphysis and upper 
part of metacarpals of proximal phalanges. The 
diameters of metaphysis, epiphysis and diaphysis of 
each EMROI were measured. The extracted features 
described the stage of skeletal development more 
objectively than visual comparison. Niemeijer et. al. 
automated the TW method to assess the skeletal age 
from a hand radiograph [19]. They employed an ASM 
segmentation method developed by Cootes and Taylor 
[20] to segment the outline of the bones. Then the mean 
image for an ROI in each TW2 stage was constructed. 
Next, an ASM was developed to determine the shape 
and location of the bones in a query ROI, so that this 
ROI can be aligned with each of the mean images in the 
third step. Then the correlation between a fixed area 
around the bones in the mean images and the query ROI 
was computed. These correlation coefficients were used 
to determine the TW2 stage in the final step. 

M. Fernandez et. al. [21] described a method for 
registering human hand radiographs for automatic BAA 
using the GP method. This method was the first step 
towards a segmentation-by-registration procedure to 
carry out a detailed shape analysis of the bones of the 
hand. Accurate results were obtained at a fairly low 
computational load. A. Fernandez et. al. proposed a 
fuzzy logic based neural architecture for BAA [22]. The 
system employed a computing with words paradigm, 
wherein the TW3 statements were directly used to build 
the computational classifier. Luis Garcia et. al. 
presented a fully automatic algorithm [23] to detect 
bone contours from hand radiographs using active 
contours. Lin et. al. proposed a novel and effective 
carpal bone image segmentation method using GVF 
model, to extract a variety of carpal bone features [24]. 
In 2005, Tristan and Arribas [25] designed an end-to- 
end system to partially automate the TW3 bone age 
assessment procedure, using a modified K-means 
adaptive clustering algorithm for segmentation, 
extracting up to 89 features and employing LDA for 
feature selection and finally estimating bone age using a 
Generalized Softmax Perceptron (GSP) NN, whose 
optimal complexity was estimated via the Posterior 
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Probability Model Selection (PPMS) algorithm. Zhang 
et. al. developed a knowledge based carpal ROI analysis 
method [26] for fully automatic carpal bone 
segmentation and feature analysis for bone age 
assessment by fuzzy classification. Thodberg et. al. 
proposed a 100% automated approach called the Bone 
Xpert method [27]. The architecture of Bone Xpert 
divided the processing into three layers: Layer A to 
reconstruct the bone borders. Layer B to compute an 
intrinsic bone age value for each bone and Layer C to 
transform the intrinsic bone age value using a relatively 
simple post-processing. Giordano et. al. 1 28] designed 
an automated system for skeletal bone age evaluation 
using DoG filtering. The bones in the EMROIs, were 
extracted using the DoG filter and enhanced using a 
novel adaptive thresholding obtained by histogram 
processing. Finally, the main features of these bones 
were extracted for TW2 evaluation. Hsieh et. al. [29] 
proposed an automatic bone age estimation system 
based on the phalanx geometric characteristics and 
carpal fuzzy information. From the phalanx ROI and 
carpal ROI, features were extracted and classified as 
phalanx bone age and carpal bone age respectively. 
Classification employed back propagation, radial basic 
function and SVM neural networks to classify phalanx 
bone age. Normalized bone age ratio of carpals was 
used to compute the fuzzy bone age. Zhao Liu and Jian 
Liu proposed an automatic BAA method with template 
matching [30] based on PSO. An edge set model was 
designed to store the middle information of image edge 
detection. The image template matching was based on 
PSO, followed by classification. TW3 classifier 
proposed by A.Fernandez et. al. (discussed in section 
3.17) was made use of to obtain the bone age. Giordano 
et. al [3 1 ] presented an automatic system for BAA using 
TW2 method by integrating two systems: the first using 
the finger bones - EM ROI and the second using the 
wrist bones - CROI. Then the TW2 stage is assigned by 
combining Gradient Vector Flow (GVF) Snakes and 
derivative difference of Gaussian filter. We have 
presented a thorough survey of literature on BAA 
methods in our previous work [32], explaining in detail 
the various work done in BAA and providing directions 
for future research. Our previous work [33] describes a 
computerized BAA method for carpal bones, by 
extracting features from the convex hull of each carpal 
bone, named as the convex hull approach. We have also 
proposed an automated BAA method to estimate bone 
age from the feature ratios extracted from carpal and 
radius bones, named as the feature ratio approach [34]. 
Our decision tree approach utilizes features from the 
radius and ulna bones and their epiphyses for BAA [35]. 
We have also exploited the epiphysis/ metaphysis 


region of interest (EM ROI) in BAA using our Hausdorff 
distance approach [36]. 

4. The Proposed system: 

The proposed system consists of two phases, 
namely: the Training phase (Fig. 3) and the Testing 
phase (Fig.4). Both the phases share the following 
modules: 

■ Image Pre-processing 

■ Morphological Feature Extraction 

■ Feature Analysis and Selection 
The last module of the training phase is the 

■ Training Module (to train the ID3 classifier) 
The last module of the testing phase is the 

■ Testing Module (to classify the image into its 
age class, thus inferring the bone age). 


Training Image 



Trained Decision Tree 


Fig. 3. Training Phase 
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Fig. 4. Testing Phase 

4.1 Image Pre-processing 

Image preprocessing is performed in two steps: 

1 . Image Enhancement 

2. PSO-Segmentation using Tetrolets 

The input image is enhanced by image smoothing 
to reduce the noise within the image or to produce a less 
pixilated image. In our system, we have done smoothing 
using a Gaussian filter to reduce noise. 

4.1.1 Edge Detection and Segmentation 

We have made use of Sobel edge detector to detect 
the edges. The Sobel edge detector uses a pair of 3 x 3 
convolution masks, one estimating the gradient in the x- 
direction (columns) and the other estimating the 
gradient in the y-direction (rows). A convolution mask 
is usually much smaller than the actual image. As a 
result, the mask is slid over the image, manipulating a 
square of pixels at a time. 

Tetrolet-based segmentation method proposed in 
this paper, makes use of tetrolets for decomposing the 
input image into sparse representation. The decomposed 
tetrolet co-efficients [38] are fed as particle solutions to 
the PSO segmentation algorithm. The algorithm 
segments the input left hand wrist radiograph and 
identifies the ROIs for further computations. 

4.1.2 Decomposition into Tetrolets: 

1. The image a ' is divide into blocks QiJ of size 4 x 4 , 

1, j = 0,...N/4 r - 1. 

2. In each block QiJ , the 1 17 admissible tetromino 
coverings c = 1, ... , 117 tire considered. 

For each tiling c, a Haar wavelet transform is applied 

to the four tetromino subsets l[ ( \ S = 0,1, 2, 3. In this 


way, for each tiling c, four low-pass coefficients and 12 
tetrolet coefficients are obtained. In QiJ, the pixel 
averages for every admissible tetromino 
configuration c = 1, ... , 117 by equation (3) and the 
three high-pass parts for / = 1,2,3 given by equation (4) 
respectively: 

a r - ,rt = [,])J =(1 witha' u, [s]= '[hmi] (3) 

=(w', , ' <c, [s|.« wilh 

(m.jOe/f* 

where the coefficients s[l,m], l,m = 0, . . . , 3, are entries 
from the Haar wavelet transform matrix: 

"1111 

J -1 -1 1 

3. The low- and high-pass coefficients of each block are 
re-arranged into a 2 x 2 block. 

4. The tetrolet coefficients (high-pass part) are stored. 

5. Step 1 to 4 is applied to the low-pass image. 

6. The tetrolet coefficients are fed as input to the PSO 
algorithm for segmentation. 

4.1.3 Overview of PSO: 

Particle Swarm Optimization (PSO) is an algorithm 
for finding optimal regions of complex search space 
through interaction of individuals in a population of 
particles. PSO algorithm, originally introduced in terms 
of social and cognitive behavior by Eberhart and 
Kennedy [39] has been proven to be a powerful 
competitor to other evolutionary algorithms such as 
genetic algorithms. PSO is a population based stochastic 
optimization technique and well adapted to the 
optimization of nonlinear functions in multidimensional 
space [40]. PSO algorithm simulates social behavior 
among individuals (particles) flying through 

multidimensional search space, each particle 

representing a single intersection of all search 
dimensions [41 1. The particles evaluate their positions 
relative to a global fitness at every iteration, and 
companion particles share memories of their best 
positions, and then use those memories to adjust their 
own velocities and positions. At each generation, the 
velocity of each particle is updated, being pulled in the 
direction of its own previous best solution (local) and 
the best of all positions (global) [42,43]. Computation of 
optimal threshold is handled here with Particle Swarm 
Optimization (PSO). There are six important control 
parameters in PSO algorithm. They are: Population 
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Size, Cognitive Learning Rate, Social Learning Rate, 
Maximum of Particle Flying Speed, Inertia Weight 
factor, and Constriction factor. The population size of 
particles refers the number of particles in iterative 
process, thus denoting components in the image here. A 
population of particles is initialized with random 
positions and velocities in d-dimensional space. A 
fitness function, f is evaluated, using the particle’s 
positional coordinates as input values. Positions and 
velocities are adjusted, and the function is evaluated 
with the new coordinates at each time-step. 

Algorithm 

Step 1: In every iteration, each particle is updated by 
following two "best" values, Personal best and Global 
best. 

Step 2: After finding the two best values, the particle 
updates its velocity and positions with following 
equation (6) and (7) 

v[ij = v[i] + cl * rand(i) * (pbest[ij - present [i]) 

+ c2 * rand( i)*(, gbest [ i ] - present / i ] ) (6) 
presen t/ij= present I ij+v[i] (7 ) 

v[i] is the particle velocity, present/ i] is the current 
particle (solution), pbest/i] and gbest/i] are defined as 
stated before, rcind(i) is a random number between (0,1) 
and cl, c2 are learning factors. Usually cl=c2=2. 

4.1.4. Implementation of PSO for tetrolet -based 
segmentation : 

The implementation of the segmentation algorithm 
consists of the following steps. 

Step 1: Swarm Formation: For a population size p, the 
particles are randomly generated between the minimum 
and the maximum limits of the threshold values. 

Step 2: Objective Function evaluation: The objective 
function values of the particles are evaluated. 

Step 3: ‘ pbest ' and ‘gbest ’ initialization: The objective 
values obtained above for the initial particles of the 
swarm are set as the initial pbest values of the particles. 
The best value among all the pbest values is identified 
as gbest. 

Step 4: Velocity computation: The new velocity for each 
particle is computed using equation (6). 

Step 5: Position computation: The new position for 
each particle is computed using equation (7). 

Step 6: Swarm Updation: The values of the objective 
function are calculated for the updated positions of the 
particles. If the new value is better than the previous 
pbest , the new value is set to pbest. Similarly, gbest 
value is also updated as the best pbest. 

Step 7: Termination: If the stopping criteria are met, the 
positions of particles represented by gbest are the 


optimal threshold values. Otherwise, the procedure is 
repeated from step 4. 

Fig. 5 provides a snapshot of the input, edge detected 
and the segmented images. 



(a) (b) (c) 

Fig. 5. (a) Pre-processed image (b) Edge detected image 
(c) Segmented image 


4.2 Feature Extraction and Analysis 

Feature extraction is performed as a way to reduce 
the dimensionality of the data. A long list of candidate 
features was calculated in order to form a powerful 
input vector. The features attempted to describe the 
morphology of the outline shape of the bones (ROIs). 
Once extracted and optimized, the vector would be used 
to train and validate the classifier. The list of candidate 
features are categorized into two groups, based on the 
ROIs they are extracted from, as: 

1 . Carpal features: (9 features) 

■ BRlength 

■ BRwidth 

■ BRdiagonal 

■ BRarea 

■ BRperimeter 

■ IN are a 

■ INperimeter 

■ BRJNratio 

■ Solidity 

2. Phalangeal features: (3 + 12x3 + 3 =42) 

■ MedianLength (3 fingers) 

■ Distal Values (3 fingers) 

o DistalLength 
o Distal Width 

o Distal Area 

o DistalPerimeter 

■ Middle Values (3 fingers) 

o MiddleLength 
o MiddleWidth 
o Middle Area 

o MiddlePerimeter 

■ Proximal Values (3 fingers) 

o ProximalLength 
o Proxima / Width 

o Proxima l A re a 

o ProximalPerimeter 
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■ DA/DP — Average of the ratios of 

DistalAreci and Distal Perimeter of the three 
fingers. 

■ MA/MP - Average of the ratios of 

Middle Area and MiddlePerimeter of the 
three fingers. 

■ PA/PP - Average of the ratios of 

Proximal Area and ProximalPerimeter of 
the three fingers. 

4.2.1 . Carpal Features: 

For the carpal bones, all the carpal bones are 
together identified as an inner region and are enclosed 
inside a bounding rectangle. Five features, namely the 
length, width, area, perimeter and the diagonal length of 
the bounding rectangle are calculated. The length and 
width are calculated as the sum of pixels along the 
boundary. 



The area, perimeter and diagonal are given by. 
Area = Length* Width 


(8) 


Perimeter = 2 x Length + 2x Width 


(9) 


Diagonal = yj Length 2 + Width 2 ^ ^ 

For the inner region of carpals, two features namely the 
inner area and perimeter are calculated as the number of 
pixels enclosed within the boundaries, and along the 
boundaries of the sample region respectively. The next 
carpal feature is the ratio between the outer and inner 
areas given by. 


BR _ IN ratio = BRarea (11) 

IN are a 

Solidity, the last carpal feature is defined as the number 
of pixels in the foreground, which quantifies the solidity 
of the carpal bones identified. Thus, a total of 9 carpal 
features are extracted. Fig. 6 provides the snapshot of 
extracting the carpal ROI features. 

4.2.2. Phalangeal Features: 

Phalangeal features are extracted from three 
fingers, namely the index, middle and third finger. For 
each finger, the three phalangeal bones: distal, middle 


and proximal are considered and from each of them four 
features are extracted (length, width, area, perimeter). 
Also the length of the median for the finger phalanges is 
calculated for each finger. The ratio of the area and 
perimeter for the distal, middle and proximal phalanges 
of the three fingers are calculated and averaged to 
obtain 3 features, DA/DP, MA/MP and PA/PP. Thus, a 
total of 42 features are extracted from the phalangeal 
bones. Fig. 7 provides the snapshot of extracting and 
cropping the phalangeal ROI. Fig. 8. (a) shows the 
extraction of features from the distal, middle and 
proximal phalanges individually and (b) shows the 
extraction of MedianLength feature. 



(a) (b) 

Fig. 7. (a) Phalangeal ROI (b) Cropped region 



(a) (b) 


Fig. 8. (a) Identifying distal, middle and proximal ROI 
(b) Extracting Length of median 

4.2.3. Feature Analysis: 

Totally 51 features are extracted from the 
segmented bones of the left hand wrist radiograph, 
producing a multi-dimensional set of data. The 
exploration of the data set can give a valuable insight of 
the expected behavior of the end system and therefore 
assist further improvement of its performance [44]. But 
the great number of features defined and computed 
makes it impossible to use them all in classification. 
Also, redundant inputs may degrade the overall 
performance of the classifier. Thus, the need for a 
feature selection procedure is very evident. In the 
proposed system. Principal Component Analysis (PCA) 
is employed for dimensionality reduction. The extracted 
features are analyzed using PCA, thus reducing the 
feature vector into 1 1 important features. PCA [45-47] 
is a vector space transformation often used to reduce 
multidimensional datasets to lower dimensions for 
analysis. Given data X consisting of N samples, PCA 
first performs data normalization by subtracting the 
mean vector m from the data. Then the covariance 
matrix Z of the normalized data (X ~m) is computed. 
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1 A 

m = —YX, ( 12 ) 

Ntf ' 

Z = (X -mxx -mf (13) 

Afterwards, the basis functions are obtained by solving 
the algebraic eigen value problem 

A = O v ' ]T<D (14) 

where 0 is the eigenvector matrix of 27, and A is the 
corresponding diagonal matrix of eigen values. Feature 
selection is then performed by keeping q (q < N) 
orthonormal eigen vectors corresponding to the first q 
largest eigen values of the covariance matrix. 

4.2.4. Feature Selection: 

The features that dominate the others in feature 
analysis are found to be the principal components. 
Those features are selected for further processing, 
discarding the lagging ones. The first seven 
predominant features selected in our system are: 
BRperimeter, Solidity, BR_INratio, MedianLength, 
DA/DP, MA/MP, and PAJPP. Thus the feature vector is 
reduced from 51 features to 7 features that contribute 
best in the bone age estimation process. The 
discrimination power of the selected features is depicted 
by plotting them against the estimated bone age, shown 
in Fig. 9. From the plots, the following are inferred: 

■ The carpal features contribute more during the 
earlier classes of the age group 1-10, while 
phalangeal features contribute much during the 
latter classes. 

■ The growth of carpal bones seems to be more rapid 
in girls, while it is gradual and slow in boys. 

■ The growth of the phalangeal bones is gradual both 
in boys and girls, at the same time their growth 
seem to be dominant in boys rather than in girls. 

■ When considering similarity, the MedianLength 
feature from the phalanges showed similar values 
for boys and girls. 

■ Much difference is encountered between boys and 
girls in the ossification of the distal phalanges, thus 
showing significant disparity in the feature DA/DP. 


4.3 Skeletal Age Inference 

This module is responsible for translating the 
features into corresponding skeletal bone age. The 
selected features after analysis cannot be directly used 
for age estimation. So they are modeled into a suitable 
format for further analysis. 

4.3.1 Feature Modeling: 

The selected features are modeled into more formal 
features to make them suitable for processing and 
classification. The features are modeled based on their 
values corresponding to age, to fall between the classes 
A to J in normal cases. If they are below the lower 
threshold bound T L , they are classified to class Z and if 
they are above the higher threshold bound T H , they are 
classified to class X. The criteria selection for the 
features is given in Table 1. Each of the selected 7 
features is modeled into any of the above mentioned 
categories, based on their values, v in years. After 
framing the formal features, classification is done using 
the renowned 1D3 classifier. 


Table 1. Criteria Selection 


S.No. 

Category 

Values v (Years) 

1. 

Z 

v <T l 

2. 

A 

T l < v < 1 year 

3. 

B 

1 year < v < 2 year 

4. 

C 

2 year < v < 3 year 

5. 

D 

3 year < v < 4 year 

6. 

E 

4 year < v < 5 year 

7. 

F 

5 year < v < 6 year 

8. 

G 

6 year < v < 7 year 

9. 

H 

7 year < v < 8 year 

10. 

I 

8 year < v < 9 year 

11. 

J 

9 year < v < T H 

12. 

X 

v >T h 
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Fig. 9. Selected Features Vs Estimated Bone age (a )BRperimeter (b)Solidity (c)BR_IN ratio (d)Median Length 
(e)DA/DP (f )MA/MP(g)PA/PP. 




Fig. 10. ID3 Decision Tree 


4.3.2. ID3 classifier: 

Studies reveal that there may be discrepancies in 
the results of classifications by different radiologists or 
at sometimes even between the classifications by the 
same radiologist in different moments [48]. This is 
because some stages are mistaken with the nearby 
stages easily. So it becomes difficult to strictly define 
the features. Instead, they have a certain degree of 
overlapping between them, since the growth pattern is 
gradual. Also the features can fall into only any one of 
the categories (i.e.) mutually exclusive. So the ID3 
classifier is more suitable for the BAA system. The root 
of the ID3 tree is formed with the feature BRperi meter , 
followed by Solidity, BR_IN ratio, MedianLength , and 
so on. We have given more importance to the carpal 
features than those from the phalanges because the 
features from the carpal bones influence more the 


estimation process for the age group 1-10. Hence carpal 
features form the root and higher level nodes of the ID3 
decision tree, followed by the phalangeal features. To 
avoid complexity, only selected layers of the ID3 tree 
are depicted in Fig. 10. The inputs to the tree will be the 
7 features modeled as values in years, and the output 
will be the skeletal class or category. The selected class 
is then mapped onto the skeletal bone age. 

4.3.3. Training Sets: 

The ID3 decision tree was trained with 6 male and 6 
female radiographs for each age group, thus with a total 
of 120 radiographs, 60 male and 60 female cases for the 
age group 1-10. After the decision tree was built and 
fine-tuned [49,50], 50 radiographs from girls and 50 
from boys were used to test it. The results showed that 
43 out of 50 were correctly classified for the girls (86%) 
and 42 out of 50 were correctly classified for the boys 
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(84%). The classification error was to assign a later 
class or a previous class. 

Results and Discussion 

The performance of the BAA system in estimating 
the bone age of patients was evaluated using a dataset 
consisting of 100 digital radiographs (50 for boys and 
50 for girls). The quality of the image is a vital factor of 
influence in the estimation process. So much 
consideration was imparted on obtaining quality 
radiographs. Also particular attention was dedicated to 
design and implement robust techniques for 
preprocessing and reliable feature extraction. During the 
preprocessing, the noise caused due to radiation and 
other external factors were eliminated. The accuracy of 
segmentation was measured using the number of under 
selected and over selected pixels. The segmentation was 
regarded as accurate if the sum of over selected and 
under selected pixels were less than 25. The 
segmentation process was accurate by 94% for boys and 
96% for girls, as tabulated in Table 1. From the 
segmented ROI, 9 carpal and 42 phalangeal features 
were extracted. The feature extraction process was 
evaluated based on the measurability, ease of 
measuring, and the reliability of the feature data such as 
area and perimeter. The accuracy of feature extraction 
from the carpal and phalangeal bones is shown in Table 
2. The ease and accuracy of feature extraction for the 


proximal phalanges and the inner carpal area were 
slightly imprecise when compared to the other features, 
the reason being that these bones were comparatively 
interior. Also in some of the radiographs, the zones 
where the phalanges are located are dark with low 
contrast. Feature space reduction was done using PCA. 
The contribution of each selected feature towards bone 
age estimation was correlated by plots (refer Fig. 11). 
Finally 1D3 classifier was used to estimate the final 
bone age. The accuracy of the classifier was measured 
in terms of perfection in classification rate. 
Classification rate, CR is defined as the ratio of the 
number of correctly classified images N C c to the total 
number of test images input to the classifier, N /P . 

N 

CR = — — (15) 

N JP 

Incorrectly classified images include the images with 
under estimated or over estimated bone age when 
compared with the bone age estimated by a radiologist 
expert. For the correct classification, a tolerance limit 
ToL of 1 year more or less was considered negligible 
(i.e. ToL = ± 1 year). The accuracy of classification is 
given in Table 3. The proposed bone age estimation 
system was accurate in estimating the bone age of boys 
by 84% and girls by 86%. 


Table 1. Accuracy of Segmentation 



No. of Images Accurate 

Inaccurate 

Accuracy 

Boys 

50 

47 

3 

94% 

Girls 

50 

48 

2 

96% 


Table 2. 

Accuracy of Feature Extraction 



BR 

(Outer) 

? rpa ! Distal 
(Inner) 

Middle 

Proximal 

Area 

100% 

86% 98% 

96% 

90% 

Perimeter 

100% 

82% 96% 

95% 

88% 


Table 3. Accuracy of Classification 



No. of Images 

Accurate 

Perfectly 

Accurate 
( < ToL) 

Inaccurate 
( > ToL) 

Accuracy 

Boys 

50 

40 

2 

8 

84% 

Girls 

50 

42 

1 

2 

86% 


Conclusions and Future work 

The work presents an efficient system for skeletal 
age assessment. The system takes digital left hand wrist 
radiographs as input and outputs the skeletal bone age 
as the output. There are two phases: the training phase 
and the testing phase. The former is used to build and 
fine tune the 1D3 decision tree classifier and the latter is 


used to estimate the bone age. The input images were 
first preprocessed by smoothening with Gaussian filter. 
Then edges of the bones were detected using Sobel 
operator and segmentation was done by a new PSO 
algorithm using Tetrolets. From the segmented ROI, 9 
carpal features and 42 phalangeal features were 
extracted. The extracted features were analyzed using 
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PC A and among them 7 dominant features were 
selected. Feature modeling was done to convert the 
selected features into a form understandable by the 
classifier. The ID3 classifier mapped the features onto 
the bone age class for the image. This bone age class 
disguised the bone age. The system was tested on a set 
of 100 radiographs (50 from girls and 50 from boys), 
achieving a success rate for bone age estimation of 86% 
for girls and 84% for boys. Future work will be focused 
on extending the system to work on the age group above 
10 years, broadening the system to include the further 
TW2 bones such as radius, ulna, etc. and integrating the 
system with PACS. 
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